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Integral equation theory of pure liquids, combined with a new "scaling approximation" based 
on a corresponding states treatment of pair correlation functions, is used to evaluate approximate 
structure factors for colloidal fluids constituted of uncharged particles with polydispersity in size 
and energy parameters. Both hard sphere and Lennard- Jones interactions are considered. For poly- 
disperse hard spheres, the scaling approximation is compared to theories utilized by small angle 
scattering experimentalists (decoupling approximation and local monodisperse approximation) and 
to the van der Waals one-fluid theory. The results are tested against predictions from analytical 
expressions, exact within the Percus-Yevick approximation. For polydisperse Lennard- Jones parti- 
cles, the scaling approximation, combined with a "modified hypernetted chain" integral equation, is 
tested against molecular dynamics data generated for the present work. Despite its simplicity, the 
scaling approximation exhibits a satisfactory performance for both potentials, and represents a con- 
siderable improvement over the above mentioned theories. Shortcomings of the proposed theory, its 
applicability to the analysis of experimental scattering data, and its possible extensions to different 
potentials are finally discussed. 

I. INTRODUCTION 

Colloidal solutions consist of mesoscopic particles, with diameters between 10 -8 m and 10 -4 m, suspended in a 
microscopic fluid. Colloidal- particles of a same chemical species are not necessarily identical, but may differ in size, 
charge or other properties au This phenomenon is called "polydispersity", and systems in which it docs not occur 
are said to be "monodisperse" . A microscopic description of a fluid with a significant polydispersity is a difficult 
task, requiring a large number of independent variables. A great simplification is, however, possible by using the 
so-called continuous-mixture, or polydisperse-mixture, formalism. Such a model, adopted also in this paper, views the 
fluid as containing an infinite number of components (p — > oo), with a continuous distribution of size and/or other 
properties. All molar fractions of the components are then replaced by a single distribution function, which describes 
the composition of the system. 

Since polydispersity can significantly affect the microscopic ordering of colloidal suspensions, it must be taken into 
account in the analysis of experimental data on such fluids. In particular, we are interested in static structure factors 
obtainable from small angle scattering of light, neutrons or X-rays. 

In this paper we present integral equation (IE) calculations for the structure of polydisperse one-species fluids 
constituted of uncharged particles with hard sphere (HS) or Lennard- Jones (LJ) interactions, with polydispersity 
in size and, for LJ systems, also in energy parameters. IEs of statistical mechanics represent a powerful, although 
approximate, tool to determine both structures and thermodynamics of fluids in a simple way. However, while using 
IEs for pure fluids or binary mixtures is a rather common and successful practice, their application to-multicomponcnt 
systems with large p or p — ► oo is problematic and, consequently, less frequent in the literature.crta Only in special 
cases, when analytical solutions are available, IE calculations for polydisperse fluids can be performed in a rigorous and 
relatively easy way. .For instance, a closed analytical formula for the scattering intensity of polydisperse hard spheres 
was obtained by Vrip from the solution of the Percus-Yevick (PY) IE. Similar expressions were derived for polydisperse 
charged hard spheres, by using the corresponding analytical solution in the mean spherical approximation!! 

Unfortunately, for most "closures" and for most potentials, including the LJ one, the IEs must be solved numerically, 
requiring more and more computer memory and time with increasing p, so that the problem soon becomes practically 
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intractable (D'Aguanno et alB reached p = 10 for Yukawa plus HS interactions). 

To overcome the impossibility of investigating polydisperse systems when IEs have to be solved numerically, alter- 
native routes have been proposed. The common idea is to replace the polydisperse fluid with an appropriate system 
with very few components. A first method, not analyzed in this paper, builds up an equivalent effective mixture with 
new components, whose molar fractions and diameters - in the case of size polydispersity - are determined 
by replacing the continuous size distribution with a p'-component histogram and requiring the equality of the first 
2p' moments of the two distributions (p' = 2 or 3 is usually sufficient). SB Other methods assume that the proper- 
ties of polydisperse fluids can be obtained, to a good approximation, from those of pure fluids, if suitably averaged 
parameters and appropriate recipes are used. This same basic idea has often been employed, in several variants 
(van der Waals equation of state, perturbation theories, etc.), to predict thermodynamics of liquid mixtures with a 
small number of components. oilil Some significant applications of these concepts have also been made to structur-al 
studies of polydisperse systems, mainly by experimentalists in the analysis of scattering data. Kotlarchyk and Chcn£2l 
introduced the decoupling approximation (DA), which is perhaps the most famous of these approaches and involves 
an exact evaluation of the form factors for all components of the polydisperse system, but approximates-every partial 
structure factor with that of an one-component fluid. A second approximation, suggested by PedersenEia and known 
as local monodisperse approximation (LMA), replaces the polydisperse fluid not with a single pure system, but with 
a superposition of non-interacting pure fluids, whose number equals that of the species in the mixture. 

The aim of the present paper is twofold. First, both DA and LMA are discussed in terms of pair correlation 
functions, to get some insight into their shortcomings. Second, we propose a simple scaling approximation (SA), 
which, with respect to DA and LMA, takes excluded volume effects more correctly into account and therefore yields 
significantly improved structure factors. By using corresponding states arguments the SA derives all pair correlation 
functions of a polydisperse mixture from an appropriate pure fluid counterpart, at the cost of only one IE computation. 
The performance of SA is tested on two typical potential models, namely polydisperse hard spheres and polydisperse 
Lennard-Jones particles. For HS systems,-|SA results for the "measurable" structure factor are compared to those 
obtained from Vrij's analytical expression,u which is exact within the PY approximation. For polydisperse LJ fluids, 
SA is tested against molecular dynamics data generated for this purpose. 



II. INTEGRAL EQUATION THEORY 



A. Basic equations 



The Ornstein-Zernike (OZ) integral-equations of the liquid state theory for p-component mixtures with spherically 
symmetric interparticle potentials area 

h a p (r) = c a p (r) + p V] Xy / dr' c ai (r') h 7 p (|r - r'|) , (1) 



"afi{r) + p^Xy / dr' c ai (r') h 7 p (|r - r'|) , 

7=1 J 



where h a p (r) = g a p (r) — 1 is the total correlation function between two particles of species a and /3 at a distance r, 
g a p (r) is the radial distribution function (RDF), c a p (r) is the direct correlation function, p = N/V the total number 
density (N — total particle number, V = volume) and x 7 the molar fraction of species 7. These equations can be 
solved only when coupled with a closure relationship, given by the exact formula 

c a p (r) = exp [-u a p (r) /k B T + j a p (r) + B a p (r)] - 1 - -f a p (r) , (2) 

plus an apurnximation to the "bridge" functions B a p (r) , which are functionals of h a p (r) and higher order correlation 
functionsoE3 (u a p (r) is the interparticle potential, ks is Boltzmann's constant and T the absolute temperature; 
j a p (r) = h a p (r) — c a p (r)). The OZ equations admit an analytical solution only in a relatively small class of cases, 
for some potentials and some peculiar closures. 1— . 

In the first case considered in this paper, i.e., for hard sphere (HS) particles with additive diametersES a a , corre- 
sponding to the potential 

/ \ — f +00, r < a a p = (a a + ap)/2 
u a p [v) - l Qj r > a ^ [6) 

an analytical solution is possible if one adds to the exact hard core condition, h a p (r) = — 1 for r < <J a p, the Percus- 
Yevick (PY) approximation 
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B a p (r) = In [1 + j a p (r)] - j a p (r) for r > a a p, (4) 

which is equivalent to c a p (r) = for r > (Tap. 

Only a numerical solution is feasible in the second case of this paper, i.e., the Lennard- Jones (LJ) potential, 



u a /3 (r) = 4 e a p 



(5) 



where (Tq,^ are LJ diameters and e a p energy parameters (well depths). The number of independent LJ parameters 
is reduced by assigning individual parameters (a a ,s a ) to each species a and obtaining the cross-interactions from 
combination rules 

a a p = (a a + op)/2 (Lorentz rule), (6) 
£q/3 = y^aSp (Berthelot rule). (7) 

The OZ equation is solved using the modified hypernetted chain (MHNC) closure described in Ref. [l5| which is one 
of the best IE approximations for both structures and thermodynamics of LJ one-component fluids. All versions of 
the MHNC theory assume that the bridge function B(r) has roughly the same functional form for all potentials and 
replace the unknown B(r) of the system under study with that of a reference system whose properties are known. 
We approximate the LJ one-component bridge function with that of an appropriate HS fluid, i.e., i?^j(r; a, p,T) ~ 
BHs{ r \d, p), where d is an equivalent HS diameter, which depends on p and T. Our_choice for Bhs is a slight 
modification of an empirical analytical approximation proposed by Malijevsky and LabiHHllla (ML), and d is selected 
by equating the second density derivative of the free energy of the LJ fluid with that of the HS reference, according 
to a prescription due to Rosenfeld and Blumo Further details can be found in the original paper O 



B. Structure factors 

The Ashcroft-LangrethS partial structure factors S a p(q) are defined as 

Sap{q) = <W + P ^/x a xp h a p(q), (8) 

where S a p is the Kronecker delta and h a p(q) the three-dimensional Fourier transform of h a p (r). Appropriate linear 
combinations of-partial structure factors define global structure factors. A first example of these is the "measurable" 
structure factorH 

p p 

Sm(q) = ^2 ^2 w a {q)wp{q) y/x a xp S a p(q), (9) 

o=l /3=1 

with 

w v {q) = Fv{q) , (10) 

where F u (q) is the scattering form factor of species v, and angular brackets, (•■■), denote, here and in the following, 
compositional averages over the distribution of particles, i.e., (Y) = 53a=i x aXa for any property Y (note that 
(w*(q)) = l). 

We assume that the scattering matter has a well-defined boundary, i.e., there is a scattering core with a well- 
defined scattering volume, not necessarily coincident with the particle volume. The former, in fact, depends on the 
particle-radiation interaction, whereas the latter is determined by the interparticle repulsions and may even not be 
well-defined, as for LJ particles (only molecules with hard body repulsions possess a well-defined volume). For spherical 
homogeneous scattering cores, the form factors are 

F (a- cr scat M cx V scatt — (g ""^*/ 2 ) (id 
Ii a{q 1 cr a )<xv a 9cr scatty 2 ' ^ ' 

where cr^ catt < a a is the diameter of the scattering core of species a, V£ catt = f (c^ catt ) 3 its volume, and ji(x) = 
(sin a; — xcosx)/x 2 is the first-order spherical Bessel function. The notation F a (q; o-^ catt ) emphasizes the dependence 
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of the form factor on <r^ catt , which in general may differ from a a . In this paper, for the sake of simplicity, cr^ catt is 
taken coincident with a a , for all species and for both HS and LJ potentials. I— . 

A second global structure factor of interest is the Bhatia-Thornton number-number structure factor obtainable 
by taking all w v = 1 in Eq. (||) 

p p 

Snn(q) = ^ ^ y/XgXp S a p(q). (12) 

a=l /3=1 

While Sm{q) represents the structure factor measured in small angle scattering experiments, Snn{o) is related to the 
fluctuations in particle numbers. 



C. Polydisperse continuous limit 

All previous formulas, written in a discrete form, refer to a finite number p of components. On the other hand, 
theoretical treatments of polydispersity with continuous distributions refer to systems with an infinite number of 
components (p — ► oo). 

For HS particles polydispersity of only one property - the diameter a - is possible. For LJ particles both a and e 
might be polydisperse,|J«n this paper, however, to simplify the LJ model, energy and size parameters will be correlated 
according to the lawE3c3 

where (a) is the average diameter, ei a \ the corresponding well depth, and z is an adjustable exponent, for which we 
take the value z — 2. Other choices will be discussed later. 

The polydisperse continuous limit of the previous discrete expressions can therefore be obtained by simple replace- 
ment rules: 

x a — > dx = /(o-)der, (14) 



$>«,..•-> / d<7/((7)..., (15) 
a •* 

where f(a)da is the probability of finding a particle with diameter in the range o~S-o~+do~, and the distribution function 
f(o~) (molar fraction density function) is normalized. Specifically, we shall use the Schulz (or gamma) distributions 

f(o) = ^TC-o-'a-O* (a > 1), (16) 
T(a) 

where T is the gamma function Jh! and the two parameters a and b can be expressed as a — 1/s 2 and b — a/ (a), in 

terms of the mean value (a) and the relative standard deviation s = \J (a 2 ) — (a) 2 / (a) . The dispersion parameter s 
measures the degree of polydispersity, and varies in the range < s < 1. For s — > 0, the Schulz distribution reduces 
to a Dirac delta function centered at (a) (monodisperse limit). For small s values f(o~) is very similar to a Gaussian 
distribution (without its drawback of unphysically negative diameters). For s closer to one, f(cr) becomes asymmetric, 
with a long tail at large diameters.O The first three moments of the Schulz distribution are: (a) , (tr 2 ) = (l + s 2 ) (<r) 

and (a 3 ) - (1 + s 2 ) (l + 2s 2 ) (a) 3 . 

Whenever analytical integration is impossible, numerical integration brings back to discrete expressions (with large 
p, of order 10 2 — 10 3 ), and therefore the replacement rule of Eq. ( |l5|) becomes unnecessary. Thus no integral is needed 
in the formulas and the discrete notation is always employed, implicitly assuming x a = f(o~ a )Ao~, which is the discrete 
analogue of Eq. ( |l~4| ) (Act is the grid size in the numerical integration). 
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III. CORRESPONDING STATES AND SCALING APPROXIMATION 



To introduce our scaling approximation for global structure factors, a corresponding statesc3'EII approach will be 
used. The correspondence principle applies to systems which have conformal pair potentials, i.e., potentials of the 
same shape. The principle takes its simplest form when the potential u a of each species a, in a set of conformal 
substances, depends on two parameters only and can be written as 

u a (r) = e a u* (j^J > ( 17 ) 

where a a and e a are a characteristic length and a characteristic energy, respectively, while u* is a dimensionless 
function of the dimensionless distance r* = r/a. Such a form of u a {r) implies that all properties of a set of conformal 
fluids can be written in terms of dimensionless reduced variables, e.g. temperature T* = kBT/s a , number density 
P* a = P a li and pressure p* a = p a al/e a . 

When written in terms of reduced distance r* = r/o~ a and wavevector = qo~ a , the RDF and structure factor of 
any pure fluid of species a can be derived by scaling as 

g a {r;p,T;cr a ,e a ) = g(r* a ;p* a ,T*), (18) 
S a (q;p,T;cr a ,e a ) = S (q a ;p* a ,T*) , (19) 



where g and S are functions common to the entire set of conformal substances. Eqs. (|T^) and ([l9]) indicate that t 
scaling correspondence applies not only to thermodynamic variables, but also to "positions" in r— and q— space. 
Both HS and LJ potentials satisfy the scaling condition given by Eq. @ (in the HS case, since e a — 1, RDFs and 
structure factors do not depend on T). 

Conformal mixtures are those in which all pair potentials are conformal to each other and to that of a pure 
(monodisperse one-species) reference fluid, according to the scaling relation 

u af} {r) = e af3 u mono (A Q/3 r mono ) , (20) 

where r mcmo = r / cr mono, while £ a p and \ a p are parameters characteristic of the pair a,/?. For mixtures conformality of 
potentials does not imply conformality of RDFs in the same simple way as for pure fluids. Nevertheless, corresponding 
states arguments_haye sometimes been employed by postulating approximate conformality relations between mixture 
and pure RDFsO'Ej The same approach is also followed in the present paper. 

We call scaling approximations (SA) those assuming approximate conformality of all RDFs of a mixture, according 
to the relation 

9aj3 y'i P> x ; T\ {c-yA"}, {e^a}) — <?mono (^a^monoi Pmono) ^mono) > (21) 

where x, {cr^s}, {s-ys} represent the complete set of molar fractions and potential parameters, 

* _ 3 rp* _ 

rmono — r^monol mono — _ 5 V ) 



while cr mono and e mono are suitably chosen average-jDotential parameters. Generalizations with cr mono and e m0 no 
replaced by pair-dependent parameters are possiblejiil but we will restrict ourselves to the simplest case. The value 
of each g a p at r is obtained from a single functional form, appropriate to a pure HS fluid, by evaluating it at a scaled 
pair-dependent distance, r' a a = X a pr, and at a corresponding thermodynamic state. 

Since the Fourier transform of h mono (X a/3 r) is A~| h m<m o(\~pq) and S moao \ = ^mono (^a^mono) ! witn 

^mono = Q'^mono; 

(23) 

the approximate partial structure factors become 

S a f3(q) — &afi + \JXaX(3 A Q | S mono ( A a ^? mo n o ; Pmono! ^monoj — 1 • (24) 
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IV. HARD SPHERE POTENTIAL 



For clarity, the various approximations examined in this paper are presented starting from a particular physical 
system, namely the HS potential. Both DA and LMA will be slightly reformulated, to point out the underlying 
approximations in terms of pair correlation functions g a p (r) . 



A. Decoupling approximation 

Kotlarchyk and Chenlll proposed the decoupling approximation in a rather general form, to treat both polydisperse 
fluids and systems of non-spherical particles. They performed two basic approximations: 

a) First, they replaced orientation-dependent interparticle potentials with spherically symmetric ones. Strictly 
speaking, this is the actual "decoupling approximation" , which allows one to break the ensemble average present in 
the exact expression of the scattering intensity into two factors, neglecting correlations between particle orientations 
and positions. The result is the Fournet-Vrij expression for the scattering intensitycl with form factors averaged over 
particle orientations. For spherical particles, this first approximation is unnecessary and hence will not be exploited 
in the present paper. 

b) Second, the partial structure factors S a p(q) were approximated in terms of the structure factor of an appropriate 
pure fluid, as 

S a p(q) = S a p + p y/x a xp h a p(q) ~ 5 a p + ^Jx a x Si [S mouo {q) - 1] , (25) 

where the definition of the effective pure fluid must be completed suitably, depending on the particular physical 
system. In terms of correlation functions, this approximation is equivalent to assuming that: i) all the RDF g a /3{r) 
have the same dependence on r, being equal to the RDF of the pure fluid, i.e., g a p{f) — ffmonoMi or 

h a p{r) ~ h mono (r); (26) 

ii) the number density of the pure fluid is equal to the total number density of the mixture, i.e., p m ono = P- With Eq. 
@ and S mono = 1 + pmono ft-mono, this implies Eq. @. 

For polydisperse HS, Kotlarchyk and Cher£3 defined the diameter of the effective HS fluid by choosing 

%o„„ = ( ff3 ) 1/3 . ( 27 ) 

ensuring that the volume fraction (or packing fraction) of the pure fluid 77 mon o = (tt/6) PmonoC mono is equal to the 
total volume fraction of the mixture, rj = (tt/Q) p(cr 3 ). 
The resulting DA structure factors are 

aXfi $mono(<7 m0 no' ri)~l , (28) 
S M (q) DA = 1 + (w(q)) 2 [S moao (q* wno ; fj) - l] , (29) 
SNN(q) DA — Smono (^mono! V) ■ (30) 

On the left hand side of these equations the dependence on the thermodynamic state and the potential parameters of 
the mixture has been omitted for simplicity. Moreover, here and in the following, we simply write i] mono — rj instead 
of /O* nono in 5 mono , since these quantities are proportional. 

Before concluding this subsection, some remarks are appropriate. The first one is that the approximation expressed 



by Eq. (26) was.|aLready proposed in the theory of liquid mixtures: it is known as random mixture, or random mixing, 
approximationJjiiJ and sometimes is also referred to as substitutional modela The weakness of this approach is evident 
from its RDF form, Eq. (|2^), because it ignores the ordering which takes place in the presence of different particle 
sizes. Therefore, the DA is expected to be a very poor approximation for moderate or even low size polydispersity. 
A second remark is that, for a given interparticle potential it mon o, it is possible to choose among several routes to 
evaluate Smono from the pure fluid OZ integral equation, by changing the "closure" . In the luckiest cases an analytical 
solution may be available, but a numerical solution is always feasible. For polydisperse HS, Kotlarchyk and Clients 
used the PY analytical solution for the monodisperse fluid, but other more accurate closures (as for instance, the 
Ballone-Pastore-Galli-GazzilloEj (BPGG) or the Rogers- Youngo (RY) approximations) could also be employed. A 
final remark is that DA may be regarded as a scaling approximation, obtained from Eq. ( pi| ) with the choice X a p = 1 

and CTmono = (ct 3 ) 1 3 . 
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B. Local Monodisperse Approximation 



This approximation was originally formulated by PederserS for polydisperse hard spheres, but can be easily ex- 
tended to polydisperse fluids with different potentials. According to the original presentation, in LMA a p-component 
mixture is approximated by a set of p non-interacting pure subsystems (a subscript "mono-a" will be used to charac- 
terize that of species a) , and the scattering intensity is calculated as a superposition of the scattering intensities from 
the subsystems, weighted according to the size distribution of the mixture. 

For HS, LMA may be expressed, in terms of partial structure factors, as 



S a p(q) - { £ / * s a ^ « (31) 



which implies that 



p 

s M (q) LMA = 5> Q <b) s mono (q* a ,v), (32) 

a=l 

SNN(q) LMA is obtained by putting all w a (q) — 1. The pure subsystem of species a consists of hard spheres with 
diameter <T mon o-a = c Q and at a number density /9 mon o-a = p( a ' 3 ) / "'ai which differs, in general, from the density 
p a = x a p of that species in the mixture. Such a choice for the p mono - a values of the p subsystems ensures that the 
volume fraction of each of these pure fluids, ?y m ono-a = (^/6)p»nono-acr^, equals the total volume fraction r\ of the 
mixture. For S , m0 ii Pedersen used the PY analytical expression.EJ 
In terms of pair correlation functions, LMA may be written as 

P y/X a Xp h a( }{r) ~ 6 a p jOmono-Q Vono (r*J Tf) , (33) 

which shows that LMA neglects all interactions between particles with different diameters, i.e., h a p(r) = if a ^ (3. 
Pedersen justified this approximation on the ground of a physical picture, in some sense complementary to DA, which 
assumes that particle sizes and positions are completely correlated. This means that particle size varies slowly with 
position, so that every particle is surrounded by particles of the same size and the system looks locally monodisperse. 
On the other hand, LMA may simply be regarded, in a corresponding states framework, as a conformality assumption 
not for g a p but for S a /3, done to reduce the double sum to a single sum in both structure factors. 



C. Scaling Approximation 

After analyzing DA and LMA in terms of pair correlation functions, it becomes evident that excluded volume effects 
are not taken into account correctly by these approximations, since the exact hard core conditions, g a p (r) = for 
f < are not satisfied. To avoid this defect and obtain reasonably accurate RDFs of HS mixtures from pure fluid 
ones with a limited effort, we propose a scaling approximation, derived from Eq. ( |2l| ) with the choice 

1/3 

Aq/3 = (Tmono/o' a /3 and cr mono — (a ) . (34) 
Since A a/ 3r^ ono = r*^, with the definition r* af3 = r/cr a p, our SA for HS can be written as 

) • (35) 

Note that ffmono is the same as in DA, while the choice for X a fj ensures that, when r < o" a p, one gets r*^ < fJ m0 no and, 
consequently, g a p{r) = 0. Since excluded volume effects are very important for the structure of condensed fluids, it 
is therefore reasonable to expect that SA is better than both DA and LMA, although it incorrectly assumes that all 
RDF values at contact are equal: g a i3{^ap) — <7mono(cmono)- Once again the choice for cr mono ensures that the RDF of 
the pure fluid is evaluated at the same packing fraction of the mixture. From Eqs. (^|), ( pi| ) and ( |34| ) one then finds 

P P a 3 
SAi(q) SA = 1 + ^^i tt ^ w a {q)wp{q) -^y S^ono^; v) ~ 1 , (36) 

with the definition g* * = qo~ a p. SNN(q) SA is again obtained from the expression for SM(q) SA by putting all w u (q) = 1. 



7 



ur SA closely resembles, albeit it is not identical to, the so-called van der Waals one-fluid 
L ^ In vdWl thermodynamics calculations the mixture is replaced by a single pure fluid, with 



It should be noted t 
(vdWl) approximation 
averaged potential parameters cr x ,e x given by van der Waals rules: 



£x = cr 



(37) 



(of course, there is no e x and no dependence on T in the RDFs of HS fluids). Although often presented in the 
literature with a different notation, vdWl may be regarded as a scaling approximation, obtainable from Eq. (|2l]) with 
A a /3 = cr mono /cr a p and cr m0 no = cr x . This choice for CT mono implies that the vdWl RDF of the pure fluid is evaluated 
at a packing fraction, n x = (ir/6)pcr 3 , which differs from that of the mixture, n = (tt/6)p (cr 3 )- As a consequence, for 
HS we can write 



p p 



Sm(q) 



vdWl 



a=10=l ° x 



S m ono{o a p]rj x ) 1 



(38) 



D. Exact PY solution 



The closed analytical expression for Sm(q) of polydisperse HS in the PY approximation, Sm(q) can be found 
in Vrij's original papeia or in Ref. || (with all charges set to zero). The corresponding expression for SNN(q) PY is 
simply obtained by putting F v {q) = 1 for all form factors in the Sm{o) py formula. 



E. Numerical results 



We tested the results for Sm(q) obtained from DA, LMA, SA and vdWl against Vrij's Sja(q) PY , which is exact 
within the PY approximation. Such a.-cprnparison requires the evaluation of S mono (q^ iono , r) mono ) PY , for which a 
simple analytical expression is availableBd 

Using the Schulz distribution to represent the size polydispersity, the packing fraction of the HS mixture may be 
written as 



n = (tt/6)p (cr 3 ) = (vr/6V* (1 + s 2 ) (1 + 2s 2 ) , (39) 

where we have chosen (cr) as the unit of length, and defined the dimensionless density p* = p (cr) 3 . 

Since working with dimensionless variables is very convenient, we also define r* = r/ (cr), ct* = a / (a) and q* = q (a). 

The effect of polydispersity on structure factors may be studied by varying s with either n or p* constant. Quali- 
tatively, the results at a fixed total density and those at a fixed packing fraction are similar. In Figures [j] and g we 
present some results for Sm(q) at fixed packing fraction 77 = 0.3, with varying polydispersity, s = 0.1, 0.3 and 0.5 
(most of the experimental s values lie in the range 0-^0.3). The three distributions of diameters were discretized with 
a grid size Act* = 0.02, and truncated at er* ut = 1.68, 3.48 and 5.90 (where /(ct)Act w 10~ 8 ). These values of cr* ut 
correspond to polydisperse mixtures with a number of components p = 85, 175 and 296, practically intractable with 
the available algorithms for solving IEs numerically. 

At a fixed 77, calculating iSm(9) or Sj\rN(q) in the DA, LMA and SA approximations requires the knowledge of 
_ 1/3 — 

Smcmo(q, v) PY B,tq = q (cr 3 ) , qcr a and qa a p- In other words, for each q value, a single evaluation of S mon o is needed 

for DA, a number p of evaluations is required for LMA, and p(p + l)/2 evaluations for SA and wdWl. For wdWl, 

Smono is computed at ?7 m ono = Vx, rather than ^ mon o = V- To save computer time in LMA, SA and vdWl calculations, 

we avoid the repeated evaluations for each q by taking advantage of the fact that Smono does not depend on q and 

cr separately, but only on their product q = qa a p = 9*0"*^. Thus we can choose a suitable grid size Aq* and a 

number of points TV, and calculate S mono (qi, rf) PY at the grid points % = iAq* (i = 0, . . . , N — 1) only once, storing 
all values in an array. Of course, the grid points % do not exhaust all the required qa a p values. Nevertheless, if the 
grid size is small enough, the value of the continuous function S mono at q = ?*ct*^ can be approximated with that at 
the nearest grid point, whose index in the array is simply determined from the ratio q* o~* a p/ Aq* . In the worst case, 
s = 0.5, to get Sm(?) in a range < q* < <z max ~ 20 (a reasonable choice), S mono (q) must be evaluated in a range 
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< q < max(q*<7* ut ) w 120. To satisfy such a condition and get a good grid size in q*-space, Aq* = q^^/N , we 
choose q* = 7r/Acr* = 507r and A/" = 4096. |— , 

Figure [j] shows DA and LMA results for Sm(q) versus the exact PY ones. Similar plots are reported by Pedersen.li2l 
With increasing polydispersity, Sm(q) increases in the low-q region, its first peak is reduced and shifted to smaller q 
values, and the subsequent oscillations are progressively washed out. The failure of DA, even at low polydispersity, is 
evident not only at small q, but also in the first peak region. LMA is significantly better at small scattering vectors, 
but does not reproduce the shape of the first peak correctly. 

The corresponding SA and vdWl results are plotted in Figure Now the agreement with the exact PY data 
is surprisingly good in both cases: the position of the first peak is well reproduced, and its height is only slightly 
underestimated. However, SA is globally superior: for s = 0.3 and 0.5 its worst discrepancies are found near the 
origin, at q* < 2, whereas the vdWl curves are somewhat shifted with respect to the exact ones on the left side of 
the first peak. It should be recalled that, unlike vdWl, SA is evaluated at the true packing fraction of the mixture. 

Finally, it is worth looking at the number- number structure factor SpfN{q)- In Figure [5] a comparison is made 
among SNjy(q) PY , SM{q) PY , and S mono (q) PY . The results for Snn(q)^ were obtained from our closed analytical 
expression, which is exact within the PY approximation (Griffith et a/H3 presented similar data, calculated without 
using a closed formulalS ). For s = 0, all these structure factors coincide, but with increasing s the differences become 
larger and larger. In particular, SNN(q) PY exhibits a more rapid flattening of the first peak than SM(q) PY , while its 
increase in the low-q 1 region is much more dramatic (see the extreme case s — 0.9, also included). Since that Snn{Q) 
measures the fluctuations in the total particle density (irrespective of the species), this behavior of SjvaK?) provides 
more physical information than that provided by Sm(q)- Figure^ then shows the SNN{q) predicted by SA, compared 
to SNjy(q) PY . The performance of SA is good in the first peak region and beyond, but the approximation fails in 
the low-q region and is unable to reproduce the number density fluctuations correctly. This discrepancy is, however, 
less important in Sm{q)i which is the structure factor more directly comparable with small angle scattering data. In 
fact, multiplying S a p(q) by the product of form factors w a (q)wp(q) reduces the mentioned defect (as already seen in 
Figure |^); moreover, a significant part of the region near q — is experimentally unaccessible. 



V. LENNARD JONES POTENTIAL 



From e a — et a \ (a a / (a)) z and the Berthelot rule, Eqs. ( [l3| ) and (0), it follows that 

£ ^ - k B T -T*{ (a) J ' (40) 

with T* = fcsT/£/ CT \. 

As for HS mixtures, in addition to a density p* = p (<r) 3 , it is convenient to define a second dimensionless variable 
4> = p (c 3 ), or r/Lj = 07r/6, which plays the same role as the HS packing fraction r\ (although, rigorously, a LJ particle 
has no definite boundary and volume). For the Schulz size distribution, <f> — p* (l + s 2 ) (l + 2s 2 ). 

We performed IE and MD calculations for LJ mixtures at fixed T* = 1, = 0.8 {rjLj — 0.42) with z — 2 and 
polydispersity parameter s = (monodisperse case), 0.1 and 0.3. For s = 0, one has <f> f~.p* , and the corresponding 
thermodynamic state, p* = 0.8, T* = 1, is near the triple point in the LJ phase diagram.Ej With increasing s at fixed 
4>, p* decreases. Taking into account the superiority of SA with respect to DA, LMA and vdWl approximation in the 
polydisperse HS case, our IE calculations for the LJ potential were restricted onlv-to-SA (some information about the 
performance of vdWJ-for LJ potentials can be found in works by Hoheisel et o/.,lSE3 but only for binary mixtures). 

Koike and GlandtotJ presented similar Monte Carlo results for a polydisperse LJ model of n-paraffins. Unfor- 
tunately, since these authors started from an activity distribution, rather than a composition distribution, their 
simulations are not directly comparable with our IE calculations. 



A. Molecular Dynamics simulations 



The MD simulations employed 1000 particles, in a cube with periodic boundary conditions, and, using Andersen's 
stochastic collision method,E3 described a fluid constrained at constant tfi = 0.8 and maintained in contact with an 
heat bath. The equations of motion were integrated using the velocity Verlet algorithm PIW The atomic diameters a 
where drawn from a Schulz distribution at the desired s, with well depths given by Eq. (40). To reduce the sampling 
noise, the distribution was periodically regenerated during the simulation. To avoid truncation effects in the Fourier 
transforms, Sm(q) was evaluated in the simulation using directly the definition 







£i-Fi(?) exp(iq- r JI 2 
*M«) = , (41) 



where q is the wavevector, q — |q|, is the position of the i— th particle and its form factor, Eq. Jll| ). A similar 
expression was used to obtain Snn{q), and we also calculated its r-space counterpart, the number-number correlation 
functioned gjvjv(r). 

B. Scaling Approximation and numerical results 

The scaling approximation for the LJ potential is derived from Eq. ( pl| ) with the choice \ a p = cr mono /a a f3 and 
CTmono = (cr 3 ) 1 ^ 3 (analogous to Eq. fl34| ) for the HS potential), and 



(42) 



which is consistent with Eqs. ([l3]) and (f40|). The SA for LJ therefore reads 

g a0 (r; p, x, T; {er 7( s}, {e 7 a}) ^ Jmono (r*^; rj LJl T^ ono ) . (43) 



Note that T^ ono decreases with increasing s: with Eqs. (|22|), ( |42D and the Schulz distribution, T* lono = 
T* I [(l + ,s 2 )(l + 2 S 2 )] 2/3 . 

To evaluate 5 , mono (q* /3 ; t]lj, Tj^ ono ), we solved numerically, for each T^ ono value, the MHNC-ML integral equation 
with Af = 4096 and Ar* = 0.02 (which corresponds to gj^ ax = ir/Ar* — 50ir, just as in the HS case). From the 
resulting S mono , the polydisperse structure factors Sm(q) SA and SjyN(q) SA were finally evaluated from the obvious 
analogue of Eq. (Bq), in the same manner presented in detail for the HS potential. 

In Figures || andg the Sm(q) SA and Snn(q) SA results are tested against the corresponding MD data. The overall 
agreement is good. For s = 0.1 the positions of maxima and minima and the shape of both structure factors are 
correctly predicted; only the height of the first peak is slightly underestimated. For s — 0.3 the first peak is well 
reproduced in Snn(q) sa , but not equally well in Sm{q) SA - On the contrary, SNiy(q) SA is completely off in the region 
near the origin, where the discrepancy in Sm(q) SA is less dramatic. Finally, Figure |?] shows the behavior of gNN(r)'- 
also in r— space SA works quite well. In the IE calculations we employed the expressionEil 

v v 

9NN(r) = ^2^2 x a xp g af 3{r). (44) 

a=l/3=l 

The MD data for the monodisperse case, s — 0, have beeruiound to be in perfect agreement with the Monte Carlo 
results for g(r) published by Llano- Restrepo and Chapmanta (not reported in Figure [?]). 

The slightly worse performance of SA in the LJ case with s = 0.3, with respect to the HS one with the same 
polydispersity, depends only on the higher packing fraction (t/lj = 0.42, whereas rj = 0.3 for hard spheres). 

Finally, it is worth mentioning that we also performed MD simulations with different z values (z = 1 and 3), and 
tested variants of the SA with alternative choices for £ mon o (for instance, e m ono = £ {<t)7 or e mon o = £{<r) (cmono/ (c)) )• 
Only very small changes of peak heights were found in all cases. This relative insensitiveness to the well depths may 
perhaps be explained by the well known fact that at high packing fractions the structure depends essentially on the 
repulsive part of potentials, while the attractive forces sensibly affect only thermodynamics. 

An improved SA could be obtained by replacing T^ ono in Eq. ( J43] ) with T*a = ksT/eap. This choice would be 
akin to the "mean density approximation" Jf3 used for binary LJ mixtures .@0El However, we have not attempted its 
application to polydisperse fluids, since it would require the evaluation of g mon o at many different reduced tempera- 
tures, which is inconvenient when IEs have to be solved numerically. In the particular case of tha-LJ potential one 
could take advantage of an explicit analytical parametrization for the RDF proposed by Goldmano This extension 
goes beyond the scope of the present paper and might be the subject of future work. 

VI. CONCLUSIONS 

In this paper we have shown that the structure of simple polydisperse fluids can be successfully predicted by integral 
equation methods, even when no analytical solution is available and a numerical one is unpracticable because of the 
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very large number of components. We have obtained rather accurate structure factors by solving IEs for appropriate 
pure fluids and applying a recipe based upon a corresponding states approach, which assumes conformality of all 
RDFs in the polydisperse mixture. Thus the hypothesis that all RDFs have essentially the same shape, but are scaled 
with respect to each other, turns out to be physically sound. 

The good performance of the scaling approximation for rather concentrated systems is due to a correct treatment of 
excluded volumes. This feature appears to be lacking in both the decoupling approximation and the local monodispcrse 
approximation, when these are reformulated in terms of pair correlation functions. The importance of evaluating the 
properties of the reference pure fluid at exactly the same packing fraction of the mixture is also to be emphasized. 
This condition is not satisfied by the very similar vdWl approximation. A shortcoming of SA is certainly the equality 
of all RDFs at contact, but this error appears to be a higher order effect, at least at high packing fractions. 

In conclusion, the proposed SA approximation offers to small angle scattering experimentalists a simple and valuable 
tool to predict structure factors or to fit data, and has proved superior to both DA and LMA. We have shown that 
the SA is applicable to different potentials for nonionic fluids. An extension of SA to polydisperse ionic mixtures will 
be presented in a forthcoming paper. 
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FIG. 1. Structure factor SM(q) of hard spheres, for different degrees of polydispersity s, at fixed packing fraction r\ = 0.3 
(here, and in all following figures, curves at different s are shifted upwards to avoid overlapping). Solid lines: Vrij's exact PY 
results. Dashed lines: local monono disperse approximation, LMA. Dotted lines: decoupling approximation, DA. 

FIG. 2. Structure factor Sm(q) of hard spheres, for different degrees of polydispersity s, at fixed packing fraction r\ = 0.3. 
Solid lines: Vrij's exact PY results. Dashed lines: scaling approximation, SA. Dotted lines: van der Waals one-fluid approxi- 
mation, vdWl. 

FIG. 3. Comparison of exact PY results for HS structure factors Sm{c[) and Snn{q), for different degrees of polydispersity s, 
at fixed packing fraction r\ = 0.3. Solid lines: exact PY results for Snn(q)- Dashed lines: exact PY results for Sm(?). Dotted 
lines: monodisperse structure factor. 



FIG. 4. Structure factor SNN(q) of hard spheres, for different degrees of polydispersity s, at fixed packing fraction r\ — 0.3. 
Solid lines: exact PY results. Dashed lines: scaling approximation, SA. 

FIG. 5. Structure factor Sm{cl) of Lennard- Jones fluids, for different degrees of polydispersity s, at fixed <f> — 0.8 and T* = 1. 
Points: molecular dynamics data. Lines: scaling approximation, SA. 

FIG. 6. Structure factor Sjvjv(g) of Lennard- Jones fluids, for different degrees of polydispersity s, at fixed <f> = 0.8 and 
T* = 1. Points: molecular dynamics data. Lines: scaling approximation, SA. 

FIG. 7. Number-number correlation functions pjvjv(r) of Lennard- Jones fluids, for different degrees of polydispersity s, at 
fixed <j> = 0.8 and T* — 1. Points: molecular dynamics data. Lines: scaling approximation, SA. 
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